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Abstract. We use the Random Dispersion Approximation (RDA) to study the Mott-Hubbard transition 
in the Hubbard model at half band filling. The RDA becomes exact for the Hubbard model in infinite 
dimensions. We implement the RDA on finite chains and employ the Lanczos exact diagonalization method 
in real space to calculate the ground-state energy, the average double occupancy, the charge gap, the 
momentum distribution, and the quasi-particle weight. We find a satisfactory agreement with perturbative 
results in the weak- and strong-coupling limits. A straightforward extrapolation of the RDA data for 
L < 14 lattice results in a continuous Mott-Hubbard transition at Uc ~ W . We discuss the significance 
of a possible signature of a coexistence region between insulating and metallic ground states in the RDA 
that would correspond to the scenario of a discontinuous Mott-Hubbard transition as found in numerical 
investigations of the Dynamical Mean-Field Theory for the Hubbard model. 

PACS. 71.10Fd Lattice fermion models (Hubbard model, etc.) - 71.27.-|-a Strongly correlated electron 
systems; heavy fermions - 71.30-|-h Metal-insulator transitions and other electronic transitions 



1 Introduction 

The Hubbard model is the minimal lattice model of spin- 
1 /2 electrons that can describe the transition from a metal 
to a correlated-electron insulator. Long before the formal 
introduction of the model, Mott [Tl used a single-band 
picture to illustrate his idea that, at half band filling and 
above some critical strength of the mutual Coulomb in- 
teraction of the electrons, the ground state must be an 
insulator with a finite gap for single-electron excitations. 

This Mott gap is readily understood in the atomic limit 
where every lattice site is singly occupied. An additional 
electron must create a double occupancy, which increases 
the energy of the ground state by a finite amount U, the 
Mott gap in the atomic limit. For itinerant electrons with 
bandwidth W and in the strong-coupling limit, U :s> W, 
adding a hole gives rise to the lower Hubbard band of 
width Whole ~ VF in the single-particle density of states. 
Adding an electron leads to the upper Hubbard band of 
width Wdoubic ~ IF, so that the Mott gap is of the order 
oi A — U — W > 0. As U /W becomes smaller, one expects 
a Mott transition at around Uc ^ W. Note that this is a 
genuine quantum phase transition [2l[3] in the sense that 
no symmetry breaking needs to occur at Uc to change the 
transport properties from insulating to metallic. 

These qualitative arguments do not, however, deter- 
mine the precise value of the critical interaction strength, 
the size of the gap, or the Landau quasi-particle weight 
as a function of the interaction strength. Even the na- 
ture of transition is not clear: are the quasi-particle weight 



and the gap continuous functions oiU/W or do they dis- 
play jump discontinuities at the Mott-Hubbard transition? 
Since the transition is fundamentally non-perturbative, 
only exact solutions can answer such questions definitely. 

Exact solutions of interacting-electron problems are 
scarce. The Hubbard model in one dimension [1] displays a 
Kosterlitz-Thouless-type transition because the Fermi-gas 
ground state is unstable against an infinitesimal Coulomb 
interaction, Uc = 0+. This nesting instability is avoided 
in the chiral or 1 //'-Hubbard model [5] , and a continuous 
transition is observed: the discontinuity of the momen- 
tum distribution goes to zero when the gap opens linearly 
above Uc = W. Thus far, the Hubbard model could not 
be solved exactly in any dimension d > 1. 

More insight into the Mott-Hubbard transition can 
be gained in the limit of high dimensions or large lat- 
tice coordination number [HUZ]- In this limit, the Hubbard 
model can be mapped onto a single-impurity problem in 
a bath of fermionic particles whose properties must be 
determined self-consistently (Dynamical Mean-Field The- 
ory, DMFT) [HIS]. The single-particle density of states for 
the corresponding single-impurity problem cannot be ob- 
tained analytically. 

An alternative is to treat the Hubbard model with a 
random dispersion relation. The Random Dispersion Ap- 
proximation (RDA) describes the paramagnetic phase of 
the Hubbard model in infinite dimensions [2lfT0] . The RDA 
to the Hubbard model cannot be carried out analytically 
either, i.e., numerical calculations on finite-size systems 
are necessary to obtain definite results. 



Satoshi Ejima et al.: Random Dispersion Approximation for the Hubbard model 



2 

Unfortunately, the DMFT and the RDA approaches to 
the Hubbard model in infinite dimensions give conflicting 
scenarios for the Mott-Hubbard transition. All numerical 
studies of the DMFT equations favor a discontinuous tran- 
sition with a jump discontinuity in the gap, whereas the 
RDA results favor a continuous transition. 

In this work, we present numerical results for the Ran- 
dom Dispersion Approximation to the Hubbard model. In 
Section [21 we introduce the Hubbard model, formulate the 
RDA approach, and give a simple example of finite-size 
effects within the RDA. In Section[3l we give the RDA re- 
sults for the ground-state energy, the single-particle gap, 
the quasi-particle weight, and other physical quantities in 
the infinite-dimensional Hubbard model. In Section 01 we 
critically discuss the two scenarios for the Mott-Hubbard 
transition and try to reconcile the confiicting results. A 
short conclusion. Section [51 closes our presentation. 

2 Random Dispersion Approximation for the 
Hubbard model 

We start our presentation with the definition of the Hub- 
bard Hamiltonian. Next, we formulate the Random Dis- 
persion Approximation (RDA) for the Hubbard model 
which becomes exact in the limit of infinite lattice co- 
ordination number. Finally, we illustrate the RDA by cal- 
culating the ground-state energy to second order in weak- 
coupling perturbation theory. 



2.1 Hubbard Hamiltonian 

We investigate spin-1/2 electrons on a lattice. Their mo- 
tion is described by 

where Cj ^ , C; ^ are creation and annihilation operators for 
electrons with spin a on site i. Here the tij are the 

electron transfer amplitudes between sites i and j, and 
ti,i = 0. 

For lattices with translational symmetry, we have ti j = 
j), and the operator for the kinetic energy is diagonal 
in momentum space, 

k;cr 

(2) 

where e(k) is the dispersion relation. The electrons are 
assumed to interact only locally, and the Hubbard inter- 
action reads 

D ^^ni,^ni^l , (3) 



where hj^cr — c] ^c^ ^ is the local density operator at site i 
for spin a. This leads to the Hubbard model [TT], 

H = f + UD . (4) 

Since we are ultimately interested in the Mott-Hubbard 
transition, we consider a half-filled band exclusively, i.e., 
a number of electrons N that equals the number of lattice 
sites L (even). Furthermore, we consider only the param- 
agnetic phase with — Ni — L/2 electrons. 

2.2 Formulation of the approximation 

2.2.1 Hubbard model with a random dispersion relation 

The density of states for non-interacting electrons is given 

by 

^(^) = iE'^(^-^W). (5) 

k 

In the limit of high lattice dimensions and for translation- 
ally invariant systems, the Hubbard model is character- 
ized by p(e) alone, i.e., higher-order correlation functions 
in momentum space factorize, e.g., |12j 

Pqi,q2(ei,e2) = -^^^(ei-e(k-|-qi))(5(e2-e(k + q2)) 

k 

= P(ei)['5qi,q2'5(ei - £2) + (1 - <5qi,q2)p(e2)] ■ 

(6) 

This is the characteristic property of the dispersion rela- 
tion in infinite dimensions [2lll0j . 

In the Random Dispersion Approximation, the dis- 
persion relation e(k) of the Hubbard model in ^ is re- 
placed by a random dispersion relation e^'^^(k), where 
each k point of the Brillouin zone is randomly assigned a 
kinetic energy from the probability distribution p{e), the 
bare density of states in ([5]) . As has been shown earlier [2j 
[To] , the RDA for the Hubbard model is equivalent to the 
exact solution of the Hubbard model in infinite lattice di- 
mensions [Bl[71[n] when long-range order is absent. 

When we choose a semi-ellipse as our bare density of 
states, the RDA for the Hubbard model becomes exact on 
a Bethe lattice with infinite coordination number T3j, 

Poi^)-^f~(^ ' (N<W2), (7) 

where W = At is the bare bandwidth. In the following, we 
shall use t = 1 as our energy unit. 

2.2.2 Implementation on finite systems 

In order to put the RDA into practice, we must calcu- 
late physical quantities for various realizations Q of the 
dispersion relation and system sizes L, i.e., we must cal- 
culate expectation values 

Xq{L)^q(ML)\X\%{L))q (8) 
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for observables X in the ground state |'fb(i))Q of Hq{L). 
Physically meaningful quantities are obtained by averag- 
ing over a large number i?^ of realizations and extrapo- 
lating to the thermodynamic limit, 



(9) 



Typically, we choose Rl = 1000 for 6 < L < 14. For dis- 
tributions with a Gaussian form, we can determine the 
average values with accuracy 0(1/L). The accessible sys- 
tem size is the limiting factor in the RDA. 

We choose anti-periodic boundary conditions for even 
L lattice sites in momentum space 

/c = y (-L-l + 2m) , TO = 1,2, ...,L, (10) 

and determine the dispersion relation e(fc) as the solution 
of the implicit equation 



w V 



2e(fc) 
W 



2arcsin( ) . (11) 



This choice guarantees that we recover po{e) from e(fc) in 
the thermodynamic limit. Next, we choose a permutation 
Qct for each spin direction a that permutes the sequence 
{1, . . . , L} into {Qo-[l], . . . , Qaii]}- This defines a realiza- 
tion of the RDA dispersion, Q = [Q|, QJ. The numerical 
task is the Lanczos diagonalization of the Hamiltonian 



(12) 



CT 1 = 1 



In this way, we obtain, e.g., the momentum distribution 
for a realization 



(13) 



which depends on momentum k only via the dispersion 
relation e(k) [BUT] (and depends on the reahzation Q). 

2.3 Second-order ground-state energy in the RDA 

As an example, we calculate the ground-state energy E = 
(H) to second order in the interaction strength numeri- 
cally. Within the RDA we have (Q = [QiQi]) 



UL 



Eq,q, {U- L) = (L)-H^ + C/2ijg)g^ (L) 



(14) 



where E^^'^ (L) is the ground-state energy of the Fermi sea, 
which is independent of the configuration Q. Eq^q^{L) is 
the leading correction, which is obtained from Rayleigh- 
Schrodinger perturbation theory as 

(2) 

\k\<Tr IpK"' 

eieQ,ik + q))0ieQ,ip-q)) 



0<q<2TT 



(15) 



where &{x) is the Heaviside step function. Here p — q and 
k + q must lie in the first Brillouin zone, i.e., p — q = 
p — q mod 27r, k + q = p + q mod 2tt, with 



fc = - (-i + 1 + 2^) 

1j 



= 0, 



p = Y + f + ; m = 0, 



L-1, 
., i-1 



q 



27m 



n = 1, 



L-1 



(16) 



In order to illustrate the importance of finite-size effects 
in the RDA, we plot histograms of the RDA data for 
Eq^Q W for = 10000 configurations and L 
20, 128, and 256 in Fig.[Tl 
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Fig. 1. Histograms of the second-order correction Sq'q (L) 
to the ground-state energy in weak coupling for (a) L = 12, 
(b) L = 20, (c) L = 128, and (d) L = 256 sites. The lines 
give the mean values of the RDA histograms, and the curves 
are fits to the Gaussian function f{E/L) — a{L) exp[—{E/L — 
h{L)f/{2c{Lf)]. 



A closer analysis of Eq. shows that the histograms 
for Eq^q^ {L) become Gaussian for large system sizes and 
that their width shrinks proportionally to 1/L. For L < 
100, the histograms are not perfectly Gaussian yet, intro- 
ducing a small systematic error in the mean value for a 
given system size. Nevertheless, the mean value can be 
determined with an accuracy that is a factor of 1/L bet- 
ter than the width of the distribution, c{L). The finite- 
size extrapolation of the second-order correction i?^^^ (L) 
is shown in Fig. [2] 

As seen from the figure, the RDA perfectly reproduces 
the result from Rayleigh-Schrodinger perturbation theory 
for the Hubbard model in Eq. (fT8|) . see below. From the 
RDA with L < 256, we find i;(2),i<256 ^ -0.0208661(9), 
compared to the analytic result E^'^^ = -0.02086614838. 
Even if we take only the results for L — 12 and L = 20 
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Fig. 2. Finite-size extrapolation of the second-order correc- 
tion S'^' (L) to the ground-state energy in weak coupling. The 
line is a second-order polynomial in 1/L. The star marks the 
analytic result in the thermodynamic limit, see Eq. (|18|l . 



into account, a fairly accurate estimate can be obtained, 
^(2).L<20 ^ _o.0208(2). 

However, outside the perturbative limit, even L = 20 is 
beyond our numerical capabilities. The RDA calculations 
which we present in the rest of this work are based on 
the exact diagonalization (ED) method which is limited 
to L < 14 lattice sites because of the memory capacity. 
Calculations for L = 16 sites using our ED code requires 
about 33 GB memory. Therefore, numerical calculations 
for Rl = 1000 configurations are prohibitively expensive. 

In principle, larger systems can be treated using the 
Density-Matrix Renormalization Group (DMRG). How- 
ever, the nonlocal nature of the hopping makes the prob- 
lem difficult in real space, and the nonlocal nature of the 
interaction in momentum space makes the problem diffi- 
cult in reciprocal space [14j . Therefore, the maximum sys- 
tem size that can be treated reliably using the DMRG is, 
at present, only marginally larger than the L — treated 
here. 



3 Physical quantities 

In this section, we present RDA data for the ground- 
state energy, the average double occupancy, the charge 
gap, the momentum distribution, and the quasi-particle 
weight. We compare our results with those from pertur- 
bation theory in weak and strong coupling and from nu- 
merical approaches to the Hubbard model in infinite di- 
mensions (Quantum Monte Carlo, the Dynamical DMRG, 
and the Numerical Renormalization Group). 



3.1 Energy and double occupancy 

The ground-state energy per site E{U)/L and the average 
double occupancy d{U) are related by 



diU) = -(D) 



1 dEjU) 
L dU 



(17) 



In weak coupling, the ground-state energy was calculated 
to fourth-order using the Brueckner-Goldstone perturba- 
tion theory [15], 



E{U) 



3n 



U 
1 



0.02086614838 



-0.0000072475 + O (U^ 



1 



d(U) = - - 0.04173229676 U 
4 

+0.00002899 U'^ + 0{U'^) 



(18) 



(19) 



In strong coupling, a diagrammatic approach based on the 
Kato-Takahashi perturbation theory provides the results 
to 11th order ^MI\, 
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Fig. 3. Histograms of the ground-state energy for U/t = 1 
(left), U/t = 4.2 (middle), and U/t = 6 (right) for L = 14 
(bandwidth W/t — 4). The lines indicate the mean values. 

In Fig. [3l we show three histograms of the ground-state 
energy distribution for U/t = 1, U/t — 4.2, and U/t = 6 
(bandwidth W/t = 4) for L — 14, calculated using ED 
with the RDA. The lines indicate the mean values. In the 
weak-coupling regime, the histograms are rather symmet- 
ric and the relative width of the distribution is quite small, 
so that the finite-size error is also small in this parameter 
regime. In strong coupling, on the other hand, the his- 
tograms are fairly asymmetric and their relative width is 
much broader than in the weak-coupling region. There- 
fore, the finite-size effects are much larger and the extrap- 
olated data are less accurate than in the weak-coupling 
regime. The inset of Fig. [4] shows the finite-size-scaling 
analysis of the ground-state energy for various interaction 
strengths. The finite-size dependence of the average values 
is weak, so that reasonable estimates for the ground-state 
energy can be obtained from the extrapolation of small 
system sizes. 
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Fig. 4. Extrapolated data for the ground-state energy as 
a function of the interaction strength in the RDA (circles), 
the DDMRG for a metallic solution (open triangles) [18], the 
DDMRG for an insulating solution (filled triangles) [19], QMC 
(stars) [171120] , weak-coupling perturbation theory [Eq. (|18p ] 
(solid line) [15], and strong-coupling perturbation theory [Eq. 
(|20p ] (dotted line) [Tfl. Inset: ground-state energy as a func- 
tion of inverse system size {L < 14) in the RDA for various 
interaction strengths. 



In Fig. [U we also compare the extrapolated RDA re- 
sults as a function of the interaction strength U/t with 
those from Quantum Monte-Carlo (QMC) [T71I20] . the Dy- 
namical Density-Matrix Renormahzation Group method 
(DDMRG) [18,19 , weak-couphng perturbation theory, see 
Eq. (|18p 15j, and strong-coupling perturbation theory, see 
Eq. (|20p 17J. As expected from the histograms and from 
the extrapolation of the ground-state energy, the RDA re- 
sults agree very well with those from all other methods 
in the weak-coupling regime. However, for U > W the 
agreement is worse, as the ground-state energy extrapo- 
lated from the RDA is consistently lower than the ground- 
state energy obtained from the DDMRG, QMC, and the 
strong-coupling expansion. 

In Fig. [5] we show the histograms for the double oc- 
cupancy for L ~ 14 sites and U/t = 1, U/t = 4.2, and 
U/t = 6 (with bandwidth W/t = 4) . As for the ground- 
state energy, the histograms are fairly narrow for small 
interaction strengths and are rather broad for large cou- 
pling. Correspondingly, the extrapolated data are more 
reliable for small than for large interaction. Interestingly, 
a double-peak structure appears for U « W. This could 
be a signature of two coexisting solutions for intermedi- 
ate interaction strengths [3] . We shall further discuss the 
significance of this structure in Section [4] 

We now extrapolate the mean value of the double occu- 
pancy as prescribed in Eq. ([9]) . We show the resulting data 
for the double occupancy in Fig. [6] where we compare our 
results with those from QMC the DDMRG [IB 

[T9] . weak-coupling perturbation theory, see Eq. ([19)) [15^, 
and strong-couphng perturbation theory, see Eq. ([2T]) [17j . 
Again, the RDA reproduces the weak-coupling expres- 
sion very well, and also provides a reasonable descrip- 
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Fig. 5. Histograms of the double occupancy in the RDA for 
interaction strengths U/t = 1 (left), U/t = 4.2 (middle), and 
U/t = 6 (right) for L = 14 (bandwidth W/t = 4). The lines 
indicate the mean values. 
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Fig. 6. Extrapolated data for the double occupancy as a 
function of the interaction strength in the RDA (circles), the 
DDMRG for a metallic solution (open triangles) [18] . the 
DDMRG for an insulating solution (filled triangles) [19], QMC 
(stars) [171120] , weak-coupling perturbation theory [Eq. (|19p ] 
(solid line) [15], and strong-coupling theory [Eq, (I2ip ] (dotted 
line) [13, Inset: double occupancy as a function of the inverse 
system size for L < 14 in the RDA for various values of the 
interaction strengths. 



tion of the strong-coupling limit. In the transition region, 
3,5 < U < 5, the RDA result for the double occupancy 
smoothly interpolates between the DDMRG and the QMC 
results for the metallic phase and the insulating phase. For 
a further discussion, see Section [4] 



3.2 Charge gap 

We calculate the single-particle gap, 

AaiU; L) = Eq{U; L + 1) + Eq{U; L - 1) - 2Eq{U- L) 

(22) 

for each configuration Q, where Eq{U ; N) is the ground- 
state energy for a system with N particles. From the 
strong-coupling theory [l6], the charge gap in the ther- 
modynamic limit is obtained as 



A(U) = U -4- - - ^ +0 ( ^ 



(23) 
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Again, the histograms for the charge gap are fairly broad. 
We mention in passing that the ground-state momenta for 
different particle numbers are not necessarily the same. 
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Fig. 7. Extrapolated data for the single-particle gap as a 
function of the interaction strength U in the RDA (circles), 
the DDMRG for an insulating solution (triangles) 18 , and 
the second-order strong-coupling expansion [Eq. (I23|l ] (solid 
line) [16]. Inset: single-particle gap as a function of the inverse 
system size (L < 14) in the RDA for various values of the 
interaction strengths. 



In Fig. [21 we compare the RDA results with the pre- 
dictions from DDMRG [TH| and the strong-couphng ex- 
pansion. Evidently, the RDA extrapolation from data for 
small systems [L < 14) overestimates the size of the gap 
and underestimates the stabihty of the metallic solution. 
Therefore, the RDA predicts a insulator-to-metal transi- 
tion at C/ w HOj. 



3.3 Momentum distribution and quasi-particle weight 

In the limit of infinite dimensions [6l[7] and within the 
RDA, the momentum distribution depends on k through 
the bare dispersion relation e(k) only, see Eq. ([T3|) . 

In Fig. [HI we show the results for the momentum dis- 
tribution for weak to intermediate couplings. For U < 
W/2 — 2, we find very good agreement between the re- 
sults from the RDA and those from weak-coupling pertur- 
bation theory to fourth order in the interaction strength. 
It is known from Feynman-Dyson perturbation theory [21j 
that the fourth-order expression becomes unreliable for 
U > 0.6W = 2.4. Therefore, we do not plot the perturba- 
tive result for U = W. 

Fig. [5| shows that the finite-size effects are strongest 
close to the Fermi energy, Ep = 0. Consequently, a reli- 
able calculation of the size of the jump discontinuity at 
Ep is difficult. Since the self-energy is independent of mo- 
mentum [6 1 I7 1 I9]. the quasi-particle weight is obtained from 
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Fig. 8. RDA data for the momentum distribution for various 
system sizes and interaction strengths U/t = 1,2,3,4 (sym- 
bols) in comparison with the results from weak-coupling per- 
turbation theory 15 (lines). 
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Fig. 9. Histograms for the quasi-particle weight for U /t = 2 
(left) and U/t = 4.2 (right) for L = 14 sites. The long vertical 
lines indicate the mean values of the RDA histograms. 



where the second line follows from particle-hole symmetry. 
For finite system sizes, we approximate the quasi-particle 
weight by the difference of the momentum distributions 
at fc = zLtt/L which are closest to e = Ep = 0, 



ZQ{U;L)^2nQ{e{-7r/L))-l 



(25) 



(24) 



and then extrapolate to the thermodynamic limit. The 
Brueckner-Goldstone J5J and Feynman-Dyson [21] per- 
turbation theories give the following result to fourth order 
in the interaction strength, 

Z{U) = 1-0.0817484 C/^+O.OOSSOlSSf/^-hO {U^) . (26) 

This expression can be used for U < 0.614^ = 2.4 ^2l!|. 

In Fig. [9] we show the histograms for the quasi-particle 
weight for U/t = 2 and U/t = 4.2 for L = 14 sites. For 
small interaction strengths, the histograms show a fairly 
narrow peak, and the mean value is well defined. For 
interaction strengths close to the metal-insulator transi- 
tion, the distribution is broad because, for finite systems. 
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even an 'insulating' configuration has a finite jump in the 
momentum distribution, Z{L;U > T^) = 0{W/{LU)). 
Therefore, we expect ZQ{L;U/t = 4.2) = 0{1/L) which 
gives Zq{L = 14:;U/t = 4) « 0.07 for our system size. 
Therefore, at intermediate to large coupling strengths, the 
finite-size effects for the quasi-particle weight are quite 
large. 
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Fig. 10. Extrapolated data for the quasi-particle weight as 
a function of the interaction strength in the RDA (circles), 
the DDMRG for a metallic solution (triangles) [15], the NRG 
(squares) |22) . and fourth-order weak-coupling perturbation 
theory [Eq. ^] [T5ll2T] (solid line). The dotted Une is a guide 
to the eye only. 



In Fig. [TO] we compare the extrapolated RDA results 
with the predictions from DDMRG NRG T2J and 
the weak-coupling expansion [15,21 . Evidently, the RDA 
extrapolation for Z{U) from the data for small systems 
{L < 14) agrees very well with all other methods for U < 
3. For U > 3, the RDA quasi-particle weight goes down 
quickly and vanishes at Uc ~ W, where the gap opens [10 . 
All other methods used for the solution of the Dynamical 
Mean-Field Theory for the Hubbard model find metallic 
solutions which exist up to [4,2 ~ 1-5W = 6. 

We shall discuss these contrasting findings in the next 
section. 



4 Discussion 



We first lay out the two conflicting scenarios for the Mott- 
Hubbard transition in the Hubbard model with infinite 
coordination number. The RDA results suggest a con- 
tinuous transition, whereas numerical treatments of the 
self-consistency equations from the Dynamical Mean-Field 
Theory support a discontinuous transition. We discuss ob- 
jections to this scenario as well as possible remedies. Fi- 
nally, we re-analyze the RDA results of Section [3] for sig- 
natures of coexistence of metallic and insulating ground 
states. 



4.1 Scenarios for the Mott transition 

Within the RDA or in the limit of large coordination num- 
ber, the metallic phase is marked by a finite quasi-particle 
weight Z{U), while a finite gap for single-particle exci- 
tations Z\(C/) characterizes the Mott-Hubbard insulator. 
Within the limits of the accuracy of the extrapolation, the 
RDA predicts the Mott transition to be continuous: as a 
function of the interaction strength {/, the quasi-particle 
weight goes to zero at Uc ~ W, where the gap opens. 
All thermodynamic quantities, such as the ground-state 
energy and the average double occupancy, are continuous 
functions of the interaction strength U. 

The scenario of a continuous transition is in conflict 
with the prediction of a discontinuous transition, as has 
been obtained from numerical solutions of the self-con- 
sistency equations in the Dynamical Mean-Field Theory 
(DMFT) [9J. A metallic solution is found for < C/ < 
Uc,2 ~ 5. 9... 6. 2, whereas an insulating solution exists 
for Uc.i w 4.5 . . .4.8 < U; the precise values for [/c,i and 
Uc,2 are stiU under discussion [I6lll7[|22[[23ll24] . Since the 
metallic solution is found to be lower in energy than the 
insulating solution, the metal-to-insulator transition oc- 
curs at Uc — Uc,2 so that the quasi-particle weight goes 
to zero continuously but the gap jumps to the pre-formed 
value associated with the insulating solution. 

A finite coexistence region at temperature T = im- 
plies a line of first-order phase transitions in the {U,T) 
phase diagram below a critical temperature Tc- It is tempt- 
ing to apply the latter to explain the phase diagram of 
V2O3 [9[|25j. which also displays a line of first-order phase 
transitions between a 'metallic' and an 'insulating' phase. 
Note, however, that the vanadium compound shows very 
large volume changes across the transition so that lat- 
tice effects probably dominate the behavior at the tran- 
sition Moreover, the critical temperature Tc derived 
from the electronic mechanism alone appears to be too 
low to account for the experimentally observed critical 
temperature. 



4.2 Problems of the scenario of a discontinuous 
transition 

The scenario of a discontinuous transition has been criti- 
cized on mathematical as well as on physical grounds. 



4.2.1 Kehrein argument 

For U U^2 7 the single-particle spectral function of the 
metallic solution consists of a quasi-particle resonance of 
width ZW ^ W which is energetically isolated from the 
lower and upper Hubbard bands, which are split by the 
pre-formed gap Ziprc- As shown by Kehrein [57], a strict 
separation of energy scales implies that the self-energy 
must display a singularity on an intermediate energy scale 
ZW <C Ec = VZW <C /\prc- Such a singular behavior 
of the self-energy is very problematic because the DMFT 
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mapping of the Hubbard model onto an effective single- 
impurity model is based on the Dyson equation and on 
the skeleton expansion around the weak-couphng limit. 
The corresponding resummations of the perturbation se- 
ries become questionable if the resulting self-energy dis- 
plays singularities as a function of frequency. 

A way out of this dilemma is the concept of a quantum 
phase transition 'as a function of frequency'. As shown by 
Metzner [28], the DMFT mapping can also be achieved 
perturbatively starting from the atomic limit. Therefore, 
one can argue that the weak-coupling perturbation the- 
ory converges for all frequencies in the region < U < 
J7c,i and that the strong-coupling perturbation theory con- 
verges for all frequencies in the region U > Uc,2- In the 
coexistence region, J7c,i < U < J7c,2j the weak-coupling 
perturbation theory converges for a finite frequency in- 
terval \uj\ < Cc, while the strong-coupling perturbation 
theory converges outside this interval. Therefore, a quan- 
tum phase transition 'as a function of frequency', signaled 
by a divergence in the self-energy at ec, is contained in 
the metallic solution. It becomes a true transition when 
Z{Uc,2) — e.c{Uc.2) — Q a.t U = Uc:2- Note, however, that 
this is a justification a posteriori^ and it is not clear a pri- 
ori that numerical methods for the solution of the DMFT 
self-consistency equation can treat such divergences prop- 
erly. 

Recently, Karski et al. [24 proposed that bound states 
exist in the coexistence regime. In their DDMRG data, 
they detected signatures of resonances at the edges of the 
lower and upper Hubbard bands. However, no such reso- 
nances are discernible in recent QMC data of Bliimer [201 ; 
see, however, Ref. j8j. While the possibility of bound states 
is not discussed in Kehrein's analysis, the existence of any 
kind of singularities in the density of states makes the 
resummation of the perturbation series disputable. 

4.2.2 Logan-Nozieres argument 

A physical argument which questions the scenario of a dis- 
continuous transition was presented by Logan and Nozi- 
eres [5S]. In the metallic phase close to J7c,2, a small num- 
ber of itinerant electrons, TVjtin — ZN, must screen a large 
number TVioc = (1 — Z)N of localized spins. For a regular 
Kondo screening process with Kondo energy Tk ~ ZW , 
the screening energy is gained once per screening electron, 
i.e., it is of the order Eg^i,, = MtinTk « Z'^WN. The itin- 
erant electrons must be excited over the performed gap 
Apre- Therefore, the loss of energy is given by i^ioss — 
A^itin^pre = ZNAp^e- Balancing the two energies thus 
gives ^prc ~ ZW, i.e., the pre-formed gap and the quasi- 
particle weight both vanish at a (single) critical Uc- 

This argument does not take into account the possibil- 
ity that the 'metallic' pre-formed gap can be larger than 
the corresponding gap in the insulating solution. Such a 
behavior has recently been observed in Ref. [23]. There- 
fore, there is an additional energy gain i^iow = ZNAp^c in 
the metallic system because the quasi-particle resonance 
pushes the filled lower Hubbard band downwards in energy 
by approximately ZApy-^. This level repulsion, rather than 



the Kondo screening, may account for the existence of 
the quasi-particle resonance in the pre-formed gap. Again, 
this is merely an a posteriori interpretation of the results 
obtained from numerical solutions of the self-consistency 
equations in the Dynamical Mean-Field Theory. 



4.3 Coexistence of metallic and insulating ground 
states in the RDA? 

Kehrein's analysis fST" points to two basic problems of 
the DMFT: (i) How many self-consistent solutions exist 
for U > 07 (ii) Which of the self-consistent solutions cor- 
responds to the Hubbard model in infinite dimensions? 

It is conceivable that multiple 'metallic' and 'insulat- 
ing' solutions of the DMFT equations exist. Some of them 
may easily be detected numerically, others may be elu- 
sive. Some solutions ought to be discarded because their 
analytical properties violate the conditions on which the 
mapping was based in the first place. Therefore, solutions 
of the DMFT equations that advocate a continuous Mott- 
Hubbard transition may be more difficult to find numeri- 
cally .23], but they could be relevant. 

It is the main advantage of the Random Dispersion Ap- 
proximation to the Hubbard model that it is not based on 
a set of self-consistency equations. Therefore, it provides 
an independent approach to the Hubbard model with infi- 
nite coordination number. Unfortunately, the system sizes 
which we can treat are fairly small, and it is difficult to 
come to definite conclusions about the existence of a con- 
tinuous or a discontinuous Mott transition. 

As mentioned at the end of Section jS] the RDA pre- 
dicts a continuous transition at C/c ~ [1^. One may 
ask if and how the RDA on finite lattices could describe a 
coexistence region. As a bare ground-state method, in the 
infinite-system limit it should give a metallic solution for 
U < Uc,2 and an insulating solution with a finite, sizable 
gap for U > Uc,2- The average double occupancy would 
follow the weak-coupling curve up to U/t = 5, e.g., its 
value at U/t = 4.2 would be d{U/t = 4.2) = 0.077. When 
we look at the RDA histogram for the double occupancy 
for L — 14: sites. Fig. (5] and at the corresponding finite- 
size extrapolation, inset of Fig. jH] we should disregard the 
data for dg < 0.053 (configurations to the left of the mini- 
mum in the histogram dQ{U)) in order to obtain the value 
d{U/t = 4.2) = 0.077 from an extrapolation of the finite- 
size data. The disregarded data would be the finite-size 
artifact of the insulating solution with its smaller values 
of the average double occupancy. 

In order to test this hypothesis, we present the his- 
tograms for the quasi-particle weight and the ground-state 
energy for U/t = 4.2 and L = 14 in Fig. [Til In the fig- 
ure, we distinguish between 'metallic' configurations with 
dQ[U/t = 4.2) > 0.053 and 'insulating' configurations 
with dQ{U/t = 4.2) < 0.053. As can be seen, indeed, 
configurations with a small value of dQ{U/t = 4.2) have a 
strong tendency to have a small value of ZQ{U/t — 4.2), 
i.e., they are rather 'insulating' than 'metallic'. Corre- 
spondingly, the quasi-particle weight will go up if we dis- 
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Fig. 11. Histograms of the ground-state energy (left) and the 
quasi-particle weight (right) for U /t = 4.2 and L = 14 sites. 
Filled bars: configurations with double occupancy dg < 0.053 
in Fig. [5l unfilled bars: configurations with double occupancy 
da > 0.053 in Fig. [5] 

regard the 'insulating solutions' in the average of Z{U] L) 
and its finite-size extrapolation. 

The 'insulating' configurations have a lower ground- 
state energy than the 'metallic' configurations. Disregard- 
ing the 'insulating' configurations would increase the RDA 
prediction for the ground-state energy and possibly im- 
prove the agreement with the QMC data. Note, however, 
that the insulating solution of the DMFT equations has a 
higher ground-state energy than the metallic solution in 
the coexistence regime. Therefore, the interpretation that 
the disregarded configurations correspond to the insulat- 
ing solution in the coexistence regime is problematic. 

In order to investigate further the meaning of 'metallic' 
and 'insulating' configurations, bigger systems with more 
configurations must be analyzed. Systems with L < 12 
do not show the clear double-peak structure for the dou- 
ble occupancy as seen for L = 14 in Fig. [5] Therefore, 
we cannot provide a separate finite-size extrapolation of 
the data for 'metallic' and 'insulating' configurations. A 
'filter' which distinguishes between 'metallic' and 'insu- 
lating' configurations would introduce some uncontrolled 
bias which would act against the construction principle on 
which the Random Dispersion Approximation is based. 



5 Conclusions 

In this work, we have used the Random Dispersion Ap- 
proximation (RDA) to investigate the Hubbard model at 
half band filling on a Bethe lattice with infinite coor- 
dination number. We have employed the numerical ex- 
act diagonalization method in real space to calculate the 
ground-state energy, the average double occupancy, the 
momentum distribution, and the quasi-particle weight for 
1000 realizations of the dispersion relation on systems with 

6 < L < 14 lattice sites for various interaction strengths. 
We reproduce earlier RDA results with a higher resolution 
because we have used twenty times more configurations 
than in previous investigations [IO l fT6 j l2T| . 



In principle, the RDA describes the Hubbard model in 
infinite dimensions adequately, but its usefulness is lim- 
ited by the small system sizes which can be treated using 
the Lanczos technique. The quantitative agreement with 
results for physical quantities in the weak-coupling and 
strong-coupling limits of the Hubbard model in infinite 
dimensions is satisfactory, but finite-size effects clearly 
reduce the significance of its prediction of a continuous 
Mott-Hubbard metal-to-insulator transition at [4 ~ VF 
{W: bandwidth). 

Numerical treatments of the Dynamical Mean-Field 
Theory for the Hubbard model in infinite dimensions in- 
dicate that the Mott-Hubbard transition is discontinu- 
ous with a finite coexistence region Uc,i < U < Uc,2 
for a metallic and an insulating solution. In the present 
study, we find that the histograms, e.g., for the double 
occupancy, are actually bimodal in some energy interval 
W < U < 1.5PF, which could be interpreted as the sig- 
nature of a coexistence region of 'insulating' and 'metal- 
lic' random configurations. In the thermodynamic limit, 
one of the two peaks should shrink to zero and reveal the 
unique metallic or insulating ground state. In order to cor- 
roborate this hypothesis, much larger system sizes must be 
analyzed. 

The RDA to the Hubbard model is free of the prob- 
lems that are generated by the mapping of the Hubbard 
model in infinite dimensions to a single-impurity problem 
whose properties must be determined self-consistently. If 
the RDA calculations could be done on larger system sizes, 
the concept of a discontinuous Mott-Hubbard transition 
could be verified independently. 
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